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Abstract 

Hypoxia is central to both ischaemic and neoplastic diseases. How- 
ever, the non-coding transcriptional response to hypoxia is largely 
uncharacterized. We undertook integrated genomic analyses of 
both non-coding and coding transcripts using massively parallel 
sequencing and interfaced this data with pan-genomic analyses of 
hypoxia-inducible factor (HIF) and RNApol2 binding in hypoxic 
cells. These analyses revealed that all classes of RNA are pro- 
foundly regulated by hypoxia and implicated HIF as a major direct 
regulator of both the non-coding and coding transcriptome, acting 
predominantly through release of pre-bound promoter-paused 
RNApol2. These findings indicate that the transcriptional response 
to hypoxia is substantially more extensive than previously 
considered. 
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Introduction 

Cells respond to physiological insults by altering their transcrip- 
tional output. Hypoxia (low oxygen level) is an important environ- 
mental stress with a central role in many physiological responses 
including adaptation to altitude, exercise, growth and development, 
as well as major pathophysiological processes such as ischaemic 



vascular disease, inflammation, wound healing and cancer [1]. The 
transcription factor, hypoxia-inducible factor (HIF) orchestrates 
many of these responses to limit the disturbance of oxygen homeo- 
stasis or promote repair processes, through the transactivation of 
protein-coding genes with key roles in pathways such as apoptosis, 
differentiation, proliferation, energy metabolism and growth factor 
production [2]. Given the fundamental role of disturbed oxygen 
homeostasis in human disease and the potential for therapeutic 
manipulation of the HIF pathway, much interest has focused on 
understanding the extent and architecture of these pathways. 

To date, pan-genomic analyses of the transcriptional response to 
hypoxia using microarrays have focussed on protein-coding RNAs 
(mRNAs) and some microRNAs (miRNAs) [3, 4]. However, recent 
genomic analyses have transformed our perspective of regulatory 
transcriptional networks. It is now recognised that mRNA forms 
only part of a much broader transcriptional output that includes 
even greater numbers of RNAs that are not translated to protein [5]. 
These 'non-coding' RNAs have important regulatory roles that fur- 
ther shape the transcriptional output of the cell. They include short 
non-coding RNA (<200 nucleotides) such as micro (miRNA), small 
nuclear/nucleolar (snRNA), piwi-interacting RNAs (piwiRNAs), and 
transfer RNAs (tRNAs) as well as long non-coding RNAs (IncRNA) 
(>200 nucleotides). Recent work has demonstrated altered non- 
coding RNA expression in many types of cancer and revealed 
functions in cell cycle regulation, apoptosis, carcinogenesis and 
metastasis [5, 6]. However, little is known about the mechanisms 
regulating the expression of the non-coding transcriptome. 

Here, we provide the first pan-genomic analysis of both non-cod- 
ing and coding transcriptional responses to hypoxia, including 
IncRNAs, miRNAs, piwiRNAs, snRNAs, and tRNAs. The analysis 
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reveals marked bidirectional hypoxia-inducible changes in transcript 
abundance across all classes of RNA and demonstrates the 
importance of HIF in regulating non-coding as well as coding tran- 
scriptional responses to hypoxia. Correlation with changes in his- 
tone marks and RNA polymerase 2 (RNApol2) positioning, indicates 
that activation of gene expression by HIF occurs commonly through 
long-range interactions that cause the release of pre-bound, 
promoter-paused RNApol2. 

Results and Discussion 

Extensive regulation of the coding and non-coding transcriptome 
by hypoxia 

To gain a global overview of the transcriptional response to 
hypoxia, we integrated expression profile and chromatin maps of 
human breast cancer MCF-7 cells incubated in either 21% ambient 
oxygen (normoxia) or 1% ambient oxygen (hypoxia) (Fig 1A). Inde- 
pendent directional and nondirectional cDNA libraries were gener- 
ated from polyadenylated RNA (polyA+) and total RNA depleted of 
ribosomal RNAs (ribo-) and analysed by high-throughput sequenc- 
ing. For transcripts present in both data sets, expression levels deter- 
mined by the two methods showed a high degree of correlation 
(supplementary Fig SI), although a subset of transcripts showed 
incomplete adenylation [7] . 

To determine transcript abundance across all classes of RNA, the 
ribosome depleted, directional RNA-seq dataset was mapped 
sequentially to databases of different RNA classes (Fig IB). 42989 
non-ribosomal transcripts were detected in normoxia and 43169 in 
hypoxia. The raw data for each class is presented in Fig 1C, whilst 
Fig ID shows the median fold-change in hypoxia, together with the 
interquartile and full ranges, after exclusion of low-abundance tran- 
scripts where this ratio might be unreliable. This demonstrates that 
all classes of RNA are regulated by hypoxia, but reveals major class- 
specific differences. Several classes of transcript (snRNAs, tRNAs 
and piwiRNAs) demonstrated overall downregulation, whereas 
mRNAs, IncRNAs and miRNAs demonstrated overall upregulation. 
To confirm these findings, representative samples of snRNA, mRNA 
and IncRNA transcripts exhibiting regulation close to the median val- 
ues were tested by quantitative (q)PCR. This analysis recapitulated 
differences between the classes, consistent with the finding of over- 
all class-specific regulation (supplementary Fig S2). In addition to 
these global changes in hypoxia, a number of transcripts in each cat- 
egory show strong up- or downregulation, ranging from 15-fold up 
to 8-fold down (piwiRNA), >1 00-fold up to 30-fold down (miRNA) 
and 4-fold up to 12-fold down (IncRNA) (supplementary Table SI). 

Several microarray-based pan-genomic analyses have described 
heterogeneous responses in the miRNA response to hypoxia 
depending on cell type, hypoxic stimulus and profiling platform - 
reviewed in McCormick et al [4]. A number of these miRNAs 
(including mir-210) were also increased in our study (supplemen- 
tary Table S2). In addition, upregulated miRNAs from our analysis 
were also induced across a panel of breast cancer cell lines (supple- 
mentary Fig S3) . 

As few IncRNAs have been included in previous pan-genomic 
analyses our data identifies IncRNAs as a new class of transcript that 
are commonly regulated by hypoxia. The majority of these tran- 
scripts remain functionally uncharacterized. However, two IncRNAs 
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Figure 1. Extensive regulation of transcript abundance by hypoxia. 

A RNA-seq was performed following 24 h culture in 21% (normoxia) or 1% 

(hypoxia) ambient oxygen. 
B Outline of the pipeline for mapping ribosomal depleted directional RNA-seq 

reads. 

C Raw data plotted as log2 fold-change by hypoxia on the vertical axis versus 
expression level on the horizontal axis. Vertical lines denote thresholds 
(higher for classes of longer RNAs) to select transcripts for further analysis. 

D Box-and-whisker plots of log2 fold-change by hypoxia for filtered 

transcripts in each class of RNA. The vertical dotted line denotes no fold- 
change. 



(NEAT1 and MALAT1) included in a previous microarray-based 
analysis [8] were the most strongly upregulated IncRNAs in our 
dataset and were also upregulated in a panel of breast cancer cell 
lines (supplementary Fig S4) . NEAT1 is a component of interchrom- 
atin-paraspeckles, which regulate nuclear retention of Adenosine to 
Inosine (A-to-I) edited RNAs [9]. MALAT1 promotes cellular prolif- 
eration, is associated with tumour growth, metastasis and recur- 
rence and is frequently upregulated or mutated in solid tumours 
[10]. In addition, H19, another hypoxia-inducible oncogenic IncRNA 
[11] was also induced in our analysis but did not have a HIF-binding 
site close by. 

Hypoxic induction of non-annotated intergenic and 
antisense transcripts 

A proportion of reads did not map to an annotated transcript of any 
class. To define whether these represent previously non-annotated 
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transcripts, they were first examined for contiguity, using the 
Cufflinks de novo transcript assembler. This identified 91 novel RNA 
transcripts with a median length of 9988 bp (range 269-38 2549 bp) 
many of which were upregulated in hypoxia (supplementary Tables 
S3 and S4) and were confirmed by qPCR (supplementary Fig S5). 

In keeping with assignment as bona fide non-annotated tran- 
scripts, reference to UCSC data (http://genome.ucsc.edu/cgi-bin/ 
hgTrackUi?hgsid = 336849981&g = cpgIslandExt) revealed that 50 
(55%) have a CpG Island within 1000 bp of the putative promoter. 
All 91 transcripts were also detected in the independently derived 
polyA+ RNA-seq dataset. Of the 91 non-annotated transcripts, 37 
did not overlap with any annotated RefSeq gene, and were therefore 
classified as non-annotated intergenic transcripts (supplementary 
Table S3). All had a low coding potential [12], indicating that they 
belong to the IncRNA class of RNA. The remaining 54 overlapped 
with a previously documented RNA, but were expressed from the 
opposite DNA strand indicating that they are non-annotated anti- 
sense transcripts (NATS) (supplementary Table S4). Anti-sense tran- 
scripts may act in cis to regulate expression of the overlapping sense 
transcript [13]. Interestingly, the majority were upregulated by 
hypoxia, induction being commonly associated either with counter- 
regulation (e.g. HIF-la, TBX2) or with co-regulation of the overlap- 
ping sense transcript (e.g. SPAG4, CELSR2) (supplementary Fig S6) . 
Thus, hypoxia-induced antisense transcripts are allied to both induc- 
tion and repression in cis adding further tiers to regulation of the 
transcriptome by hypoxia. 

Transcriptional regulation of long non-coding RNAs by HIF 

As HIF is a major transcriptional regulator of the mRNA response to 
hypoxia [14], we next examined the role of HIF in the regulation of 
other classes of RNA. As a first step we examined for spatial associa- 
tion between each class of transcript and HIF-binding sites, using 
well validated, previously published HIF-1- and HIF-2-binding sites 
in MCF-7 cells [14]. That work considered only protein coding genes 
and reported a non-random distribution of HIF binding across the 
genome. Although HIF-binding sites were strongly enriched at the 
promoters of coding genes, many HIF-binding sites were observed 
to be remote from these promoters. In the light of extensive regula- 
tion of non-coding RNAs by hypoxia, we reanalysed the distribution 
of HIF-binding in relation to the location of all transcripts. The 
number of HIF-binding transcripts (defined as the closest promoter 
to each HIF-binding site) in each category broadly mirrored the total 
number in each class. Although the numbers in many classes were 
too small to permit firm conclusions, no enrichment of HIF binding 
was observed in the vicinity of classes of RNA showing global 
downregulation (piwiRNAs, snRNAs and tRNAs) suggesting that 
they are unlikely to be directly regulated by HIF. In contrast, HIF 
binding was enriched in the vicinity of both mRNAs and IncRNAs 
(Fig 2A) with HIF-2 binding a slightly higher proportion of IncRNAs 
than HIF-1. Nevertheless, despite promoter enrichment, the majority 
of HIF-1- and HIF-2-binding sites (50% and 70% respectively) were 
still found to lie in excess of 2.5-kb from the promoter of an 
expressed transcript of any class (Fig 2B and C) . 

To assess the functional role of HIF in the regulation of different 
classes of RNA by hypoxia, we next examined for association of HIF 
binding with regulation by hypoxia across all classes of RNA using 
gene set enrichment analysis (GSEA) [15]. All expressed transcripts 
were pooled and ranked according to their hypoxic induction. HIF 



binding was strongly enriched amongst transcripts upregulated (but 
not downregulated) by hypoxia (Fig 2D). Only mRNA, IncRNA, 
miRNA and non-annotated transcripts were present in the core- 
enrichment group of HIF-binding hypoxia-inducible transcripts. 
GSEA analysis by subgroup confirmed enrichment for mRNAs 
(Fig 2E) and IncRNAs (Fig 2F), but the number of miRNAs and 
non-annotated RNAs was too small to permit statistical analysis. 
Taken together, this establishes a non-random spatial association 
between HIF-binding sites and hypoxia-inducible transcripts for 
both mRNAs and IncRNAs, pointing to a major new role for HIF in 
the transcriptional regulation of IncRNAs. 

To pursue this further, we undertook polyA+ RNA-seq analysis 
of mRNA and IncRNA in hypoxic MCF-7 cells transfected with 
siRNAs targeting HIF-la, HIF-2cx or both isoforms simultaneously. 
Approximately 25% of all IncRNA transcripts detected in the ribo- 
somal-depleted RNA-seq were detected in this analysis. Fig 3A 
portrays the results of HIF-a siRNA interventions on all HIF-binding 
IncRNA transcripts detected in this dataset. As expected, the heat- 
map confirms that many of these RNAs were upregulated by 
hypoxia and reveals a strong correlation with down regulation by 
siRNAs targeting one or both HIF-a isoforms, a result that was 
confirmed by GSEA (supplementary Fig S7). Individual datasets are 
illustrated for MALAT1 (downregulated by both HIF-la and HIF-2a 
siRNA, Fig 3B) and NEAT1 (downregulated by HIF-2o( but not by 
HIF-la siRNA, Fig 3C). Thus, the presence of HIF binding, hypoxic 
upregulation and subsequent downregulation by the suppression of 
HIF-a defines IncRNAs as a new group of HIF target genes. 

HIF-dependent recruitment and release of promoter-paused 
RNApol2 in hypoxia 

To pursue the mechanism of transcriptional activation by hypoxia 
further, we next examined RNApol2 binding and histone H3K4me3 
modification in normoxia and hypoxia. As expected, the promoter- 
associated abundance of each mark correlated with the basal tran- 
script level. Strong peaks of RNApol2 and H3K4me3 signal were 
also seen at the putative transcriptional start sites (TSS) of the previ- 
ously non-annotated transcripts consistent with assignment as bona 
fide genes (supplementary Fig S8) [16, 17]. 

Consideration of changes in RNApol2 binding in hypoxia 
revealed that for the majority of hypoxia-inducible genes, RNApol2 
was already bound at the promoter in normoxic cells. Although 
there were exceptions, at the majority of promoter sites hypoxia 
did not increase promoter-bound RNApol2 substantially. This is 
illustrated in Fig 4A which depicts changes in RNApol2 ChIP signal 
in the vicinity of the promoters of the hundred most strongly 
hypoxia upregulated transcripts of any class. In hypoxia, little 
increase in loading of RNApol2 at the transcriptional start was seen; 
instead an increase in RNApol2 is observed across the body of the 
gene. Similar patterns of change in hypoxia were seen for both 
mRNAs and IncRNAs (supplementary Fig S9). Actively transcribing 
RNApol2 can recruit histone methyl transferase activity that trime- 
thylates H3K4 either directly [18] or through p300 [19]. The distri- 
bution of H3K4me3 showed a very similar response to that of 
RNApol2, with little change in signal at the promoter in hypoxia, 
but an increase downstream of the transcriptional start site, consis- 
tent with the hypothesis that RNApol2 recruits histone methyl 
transferases that trimethylate H3K4 (Fig 4B). In contrast, consider- 
ation of DNAse hypersensitivity analyses in normoxic and hypoxic 
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Figure 2. HIF binding upregulates both the coding and non-coding transcriptome. 

A-C she proportion of HIF-la and HIF-2ci binding sites that are closest to transcribed loci of each RNA class (A). The distribution of (B) HIF-la and (C) HIF-2a binding 
sites around the transcriptional start site of the nearest expressed gene irrespective of class (black bars). The grey bars show the distribution when only active 
mRNA genes are used (analogous to previous analyses). 

D-F GSEA analysis against fold-regulation by hypoxia for (D) all HIF-binding transcripts, (E) mRNA and (F) IncRNA. 

C Log2 fold-change of transcript abundance in hypoxia, ranked in order from highly upregulated to downregulated transcripts, from which the GSEA were derived. 



MCF7 cells reported by the Encyclopedia of DNA Elements 
(ENCODE) consortium [17] revealed the expected hypersensitivity 
signals at these promoter regions, but did not reveal changes in 
hypoxia (Fig 4C). 

On close inspection the maximum RNApol2 signal was observed 
in the promoter-proximal region, 40-45 bp downstream of the tran- 
scriptional start. This suggests that promoter-proximal pausing is an 
important rate-limiting step for transcription of hypoxia-inducible 
genes and that reversal of pausing is a key step through which 
hypoxia regulates transcription [20]. To determine whether this is a 
feature of direct transcriptional regulation by HIF, as opposed to indi- 
rect effects of HIF or hypoxia itself on transcription, we performed 
similar analysis of RNApol2 distribution for genes that were down- 
regulated by hypoxia and which (consistent with indirect regulation) 
do not show promoter enrichment of HIF binding. We also compared 
hypoxia-inducible changes in patterns of RNApol2 binding and H3K4 



trimethylation amongst groups of genes that were matched for the 
extent of upregulation by hypoxia, but were selected as being either 
the closest gene to a HIF-binding site (likely to represent direct tran- 
scriptional regulation) or were located >1 Mb from the nearest 
HIF-binding site (likely to represent indirect actions of HIF or 
hypoxia). Whilst a clear change in the ratio of promoter paused 
RNApol2 to that distributed along the gene body (travelling ratio) 
was observed for genes predicted to be direct HIF transcriptional tar- 
gets this was not observed for genes predicted to be upregulated or 
downregulated indirectly by HIF or by hypoxia itself (supplementary 
Fig S10). Finally, ChlP-qPCR at a selection of these genes confirmed 
that the hypoxic induction of RNApol2 within the body of the gene 
was dependent on the presence of HIF-oc subunits (Fig 4D-F) . 

Amongst genes predicted to be direct transcriptional targets of 
HIF we observed a clear change in the distribution of travelling 
ratios in hypoxia, reflecting release of promoter-paused RNApol2 
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Figure 3. Direct transcriptional regulation oflncRNAs by HIF. 

A Heat map showing fold regulation by hypoxia and by the indicated HIF siRNA for IncRNAs adjacent to HIF-binding site that are detected in the polyA RNA-seq 
analyses. 

B, C RNA-seq and HIF ChlP-seq genome browser tracks for the two most hypoxically upregulated HIF-binding IncRNAs: (B) MALAT1 and (C) NEAT1. 



at most loci (Supplementary Fig S9E, F and Fig 4G). However, 
we also observed infrequent HIF target genes where RNApol2 
was recruited de novo in hypoxia (Fig 4H). No association was 
detected between this type of behaviour and biological function. 
However, the promoters of these transcripts also show constitu- 
tive DNAsel HSS signal indicating that other factors must limit 
normoxic RNApol2 binding at these sites. As DNA-seq cannot 
distinguish individual nucleosome positions, this might include 
mechanisms such as nucleosome positioning as well as DNA/his- 
tone methylation [21]. 

Taken together, our findings indicate that the regulation of gene 
expression by hypoxia is substantially more complex than has previ- 
ously been considered. All RNA classes were regulated bi-direction- 
ally, our analysis identifying both general responses characteristic of 
particular types of RNA and high-amplitude effects on individual 



species. Both the non-random spatial association between HIF-bind- 
ing sites and hypoxia-inducible genes and direct functional interven- 
tion by siRNA implicated HIF in the regulation of extensive 
networks of non-coding RNAs. RNApol2 profiling indicates that this 
direct regulation by HIF is mediated through activation of pre-bound 
RNApol2. These findings confirm and extend the role of HIF in the 
regulation of miRNAs and identify both intergenic and antisense 
IncRNAs as a new target class for the HIF transcriptional response. 
The ability of non-coding RNAs to regulate gene expression both in 
trans and in cis suggests the hitherto unrecognised involvement of 
HIF in a far broader network of gene regulation. Given the funda- 
mental role of hypoxia in human pathology and emerging data on 
the importance of non-coding RNAs in cancer and other diseases, 
these findings open new avenues to the better understanding of 
human disease processes. 
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Figure 4. HIF predominantly regulates transcription through release of promoter-paused RNA-pol2. 

A-C Mean distribution of normoxic (red) and hypoxic (blue) RNApol2 binding at the 100 transcripts most upregulated by hypoxia. Inset shows expanded view at the TSS 
(FPKM = fragments per kilobase per million reads). The same plots for (B) H3K4me3 and (C) DNAsel hypersensitivity are shown. 

D,E ChlP-qPCR analysis of RNApol2 within the body of (D) ALDOA and (E) NDRG1 gene shows suppression of the hypoxic induction by HIF-la & HIF-2a siRNA 
(biological duplicates). 

F Immunoblot analysis of HIF-a levels. 

C,H RNA-seq and ChlP-seq tracks illustrating (G) release of promoter-paused RNApol2 and (H) de nouo recruitment of RNApol2, despite constitutive DNAsel 
hypersensitivity. 



Materials and Methods 

Cell culture & HIF siRNA treatment 

Human MCF-7 breast carcinoma cells from the American Type Cul- 
ture Collection (ATCC) were incubated for 24 h in an In vivo2 
Hypoxia Work Station (Ruskinn Technology Ltd, Bridgend, UK) in 
an atmosphere of either normoxia (21% oxygen) or hypoxia (1% 
oxygen). Expression of HIF-la and/or HIF-2a subunits was 
suppressed as previously described [22] . 

Ribosomal depleted RNA-seq 

Total RNA was prepared using the mirVana miRNA Isolation Kit 
(Ambion; Life Technologies Ltd, Paisley, UK) and treated with 
DNasel (TURBO DNA-free™, Ambion). Ribosomal RNAs were 
depleted using RiboMinus™ Locked Nucleic Acid probes (Invitro- 
gen; Life Technologies Ltd, Paisley, UK). The resulting RNA was 
fragmented and end-repaired using the Illumina directional 
protocol. The final cDNA library was purified using Agencourt 



AMPure beads (Beckman Coulter UK Ltd, High Wycombe, UK) 
optimised to retain small RNAs. 

PolyA+ selected RNA-seq 

Directional PolyA+ RNA libraries were prepared using the ScriptSeq™ 
v2 RNA-Seq kit (Epicentre, Madison, WI, USA) and non-directional 
libraries using the TruSeqRNA Sample Prep Kit (Illumina, San Diego, 
CA, USA). 

ChlP-seq 

Chromatin immunoprecipitation was performed as previously 
described [14]. In addition to HIF ChlPs, antibodies to H3K4me3 
(Cell Signalling, Danvers, MA, USA, #9751) and RNApol2 (Santa 
Cruz Biotechnology, Dallas, TX, USA, #sc-899) were used. Libraries 
were prepared using the Illumina ChlP-seq kit. 

High throughput sequencing 

All libraries were sequenced using the GAII or HiSeq platforms (Illu- 
mina) according to the manufacturers protocol. Raw and mapped 
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data for RNA-seq and ChlP-seq is available at EMBL-EBI Array 
Express (E-MTAB-1994, E-MTAB-1995). 

Bioinformatic analysis 

After trimming adapter sequences using the FASTX-tool kit (http:// 
hannonlab.cshl.edu/fastx_toolkit/), RNA-seq reads were sequen- 
tially mapped to RNA databases and expression analysis performed 
using the edgeR (Bioconductor, Seattle, WA, USA) package. Remain- 
ing unmapped reads were examined for non-annotated reads using 
the Cufflinks de novo transcript assembler (http://cufflinks.cbcb. 
umd.edu). ChlP-seq databases were analysed as previously 
described [14]. 

Statistical analysis 

All assays were performed in duplicate or triplicate. Statistical 
analyses using two-tailed t-tests were performed in SPSS (IBM Corp, 
New York, NY, USA) and R (http://www.R-project.org). 

Supplementary information for this article is available online: 
http://embor.embopress.org 
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